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ABSTRACT 

Cosmological hydrodynamical simulations of primordial star formation suggest that the gas within 
the first star- forming halos is turbulent. This has strong implications on the subsequent evolution, 
in particular on the generation of magnetic fields. Using high-resolution numerical simulations, we 
show that in the presence of turbulence, weak seed magnetic fields are exponentially amplified by 
the small-scale dynamo during the formation of the first stars. We conclude that strong magnetic 
fields are generated during the birth of the first stars in the universe, potentially modifying the mass 
distribution of these stars and influencing the subsequent cosmic evolution. We find that the presence 
of the small-scale turbulent dynamo can only be identified in numerical simulations in which the 
turbulent motions in the central core are resolved with at least 32 grid cells. 
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1. INTRODUCTION 

Magnetic fields are ubiquitous in the local universe 
(|Beck et alj 119961) and there is gr owing evidence of 
their presence also at high redshifts (jBernet et alJl2008t 
iRobishaw et~alll2008t lMurphyH2009T ) . The seeds for these 
fields are probably a relic from the early universe, pos- 
sibly aris ing during inflation or some other phase transi- 
tion (see. lGrasso fc Rubinsteir3l2001t lSubramanianll2010l 
for reviews). Alternatively, they could be generated 
by astrophysical mechanisms lik e the Biermann battery 
(|BiermannlTl950t IXu et~aTJl2008l) or the Weibel instabi l- 
ity (jSchlickeiser fc Shuklal [200l iMedvedev et al.ll2004f >. 
Regardless of their physical origin, most models pre- 
dict weak field strengths or have very large uncertain- 
ties. Thus, magnetic fields are often ignored in studies 
of primordial star formation. 

However, recent cosmological hy drodynamica l sim- 
ulations of first star formation (lAbel et al.1 120021: 
O'Shea fc Normanll2007t lYoshida et al.l l2008: Tur k et al.1 



2001 " show the presence of turbulence in the minihalos 



where it plays an important role in regulating the trans- 
port of angular momentum. This opens up the possibility 
of exponentially amplifying initially weak seed magnetic 
fields by the action of the small-scale dynamo. Field am- 
plification via the small-scale dynamo requires both tur- 
bulent gas motions and high magneti c Reynolds numbers 
(|Brandenburg fc Subramanianll2005[ h 

The properties of the small-scale dynamo have 
been explored both in numerical simulations of 
driven t urbulence without s e lf-gravity and in analytic 
models (lHaugen et alj 12004 ISchekochihin et all 12004 



iBrandenburg fc Subramanianl I2005T ). as well as in cos 
mological simulations of magnetic field generation in 
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galaxy clusters (|Xu et al.ll2009l) . The importance of the 
small-scale dynamo during galaxy formation was previ- 
ously sugges ted by iBeck et alj (119941). lArshakian et al.1 
(|2009f) , and Ide Souza fc Qpherl (|20iOI) . It's relevance 
for primordial star form ation was previously assessed by 
ISchleicher et al.l ((2010). Possible amplification mecha- 
nisms in the protostellar disk have been previously iden - 
tifie d bvlPudritz fc Silkl (fl989l) , iTan fc Blackmanl (|200l) 
and lSilk fc Langerl (|20od L~ 



In this Letter, we show, using high resolution numer- 
ical simulations that strong and dynamically important 
magnetic fields are generated from initially weak seed 
magnetic fields by the small-scale turbulent dynamo dur- 
ing the collapse of primordial halos. The study presented 
here is a proof- of- concept investigation to show that the 
small-scale dynamo is operational during the collapse of 
gas clouds, which we apply to the conditions of first 
star formation. In Section [2] we describe our numeri- 
cal method and the initial conditions. The results ob- 
tained from our numerical simulations are presented in 
Section [3J Finally, we summarize and discuss our find- 
ings in Section |4] 

2. NUMERICAL METHOD AND INITIAL 
CONDITIONS 

We focus on the gravitational collapse and magnetic 
field amplification of the baryon-dominated inner parts 
of a contracting primordial halo, using a simplified ini- 
tial setup (rather than a full cosmological simulation) 
where we neglect the detailed thermodynamics and non- 
ideal MHD effects. The numerical simulations presented 
here were performed with the publicly avail able adaptive- 
mesh refinement (AMR) code, FLASH2.5 (|Frvxell et al.1 
2000). We solve the equations of ideal MHD including 
self-gravity with a refinement criterion guaranteeing that 
the Jeans length, 



1/2 



(1) 



with the sound speed c s and the gravitational constant 
67, is always resolved with a user-defined number of cells. 
We use an MHD Riemann solver developed for FLASH 
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that preserves positive states and proved to be highly ef- 
ficient and accurate for modeling astroph ysical problems 
involving turbulence and strong shocks (IBouchut et aT 



12 



2010) 



20071: IWaagan|[200l IBouchut et aIll2010t1Waagan et al 



The initial conditions for our numerical simulations 
were motivated fr o m la rge r-scale cosmo lo gical mod- 
els of lAbel et all d2002ti . iBromm et all <|2002t) and 
lYoshida et all (|2008l )" We initialize a super-critical 
Bonnor-Ebert (BE) sphere with a core density of pbe — 
4.68 x 10~ 20 gcm~ 3 (n B E = 10 4 cm~ 3 ) at a tempera- 
ture of T = 300 K. We also include a small amount of 
rotation of about 4% of the gravitational energy. The 
radius of the BE-sphcre is 1.5 pc which corresponds to 
a dimensionless ra dius of £ = 8.28. Note th at the criti- 
cal value is 6.451 (|Ebertlll955l: lBonnorlll956l ). The total 
computational domain is (3.9 pc) 3 in size. The initial 
conditions are furthermore characterized by the presence 
of a random initial velocity field with transonic velocity 
dispersion of amplitude 1.1 km s -1 (equal to the initial 
sound speed) and a weak random magnetic field with 
B rms ~ 1 nG. This initial field strength corresponds 
to (3 » 10 10 , where /3 is the ratio of the thermal pres- 
sure to magnetic pressure. Both, the turbulent energy 
and magnetic field spectra were initialized with the same 
power law dependence in wave number space, oc fc -2 , 
with most power on scales ~ 0.8 pc, which roughly cor- 
responds to the initial Jeans length of the c ore. Con- 
sistent with previous w orks of lOmukai et al.l (|2005h and 
iGlover fc Savinl (|2009l ). which follow the thermodynam- 
ics during the collapse, we adopt an effective equation 
of state with T = d\ogT/d\ogp +1 = 1.1 for num- 
ber densities in the range n = 10 5 — 10 10 cm~ 3 . These 
initial conditions reflect the physical parameters of the 
baryo n-dominated regim e in the centers of the first mini- 
halos (|Abel et al.ll2002T ) . Hence we can safely neglect the 
influence of dark matter at this point 5 . The applica- 
bility of the ideal MHD approximation fo r the descrip- 
tion o f the dynamics has been studied by iMaki fc Susal 
(2004) using one-zone models that follow the chemistry, 
and in particular the abundances of ionized species dur- 
ing the protostellar collapse phase. Their study suggests 
that the ionization degree is sufficiently high to ensure 
a strong coupling between ions and neutrals to maintain 
flux-freezing. 

We note that the efficiency of the dynamo process de- 
pends on the Reynolds number and is thus related to how 
well the turbulent motio ns are resolved (jHaugen et all 
120041 ; iBalsara et al"1l2004D . Higher resolution yields larger 
field amplification. To demonstrate this effect, we per- 
form five numerical simulations where we resolve Aj by 
8, 16, 32, 64 and 128 cells. 

3. RESULTS 

We find that the dynamical evolution of the system is 
characterized by two distinct phases. First, as the initial 
turbulent velocity field decays the system exhibits weak 
oscillatory behavior and contracts only slowly. Soon, 
however, the runaway collapse sets in. Fig. [1] shows a 

5 Dark matter would only have an indirect influence on the mag- 
netic field dynamics through its gravitational interaction with the 
plasma. At low densities, it may accelerate the collapse and in- 
troduce additional turbulence in the gas, which could potentially 
enhance the dynamo. 
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Fig. 1. — Two dimensional slices through the center of the col- 
lapsing core at the time when the initial field strength has increased 
by a factor of ~ 10 6 , showing the central region of our highest- 
resolution simulation (Aj resolved by 128 cells). The circle indi- 
cates the control volume Vj centered on the position of the current 
density peak. The top image shows the density and the velocity 
component in the xy-plane, indicating radial infall in the outer re- 
gions and turbulent motions in the inner core. The bottom image 
depicts the total magnetic field amplitude and the local magnetic 
field direction. 



snapshot of the central region of the collapsing core in 
our highest-resolution simulation at a time when the cen- 
tral density has increased by a factor of ~ 10 6 . The mag- 
netic field strength has grown by a factor of 10 6 , reaching 
a peak value of about 1 mG at r ~ 12. The top image 
shows the density and the velocity structure, and the 
bottom image shows the magnetic field strength and the 
local magnetic field directions. 

To understand the behavior of the system more quan- 
titatively, we need to follow its dynamical contraction. 
First, we note that the physical time scale becomes pro- 
gressively shorter during collapse. We therefore define a 
dimensionless time coordinate r, 



dt/t s {t) 



(2) 
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Fig. 2. — Time evolution of the radial density profile (panel 
a), the radial profile of the mean value of the infall velocity, v T 
(panel b), the radial profile of the rms value of v T together with 
the rms fluctuations of the lateral velocity v\, i.e., the component 
perpendicular to v T (panel c), the radial profile of the rms value of 
the total velocity (panel d) and the Mach number Ma (panel e). 
The vertical line indicates the Jeans radius at r = 12. The data 
are taken from our highest-resolution run. 

which is normalized in terms of the local free-fall time, 

1/2 



tff(t) = 



3tt 



32 GfHp.it) 



(3) 



where Pm{t) is the mean density of the contracting central 
region. We also define a control volume for gravitational 
collapse based on the Jeans volume, Vj = 47r(Aj/2) 3 /3. 
We obtain all dynamical quantities of interest as averages 
within the contracting Jeans volume, which is centered 
on the position of the maximum density. This approach 
ensures that we always average over the relevant volume 



for collapse and field amplification. 

The density profile always maintains a flat inner re- 
gion within the size of the local Jeans length. This is 
illustrated in Fig. [2^ at three different times, where a 
flat core of the size of Aj is sustained throughout the 
collapse. The radial infall motions dominate the total 
velocity in the envelope as shown in Fig. [2Jd. The tur- 
bulence is maintained on scales below Aj and dominates 
the dynamics inside the core region as can be seen from 
Figs.^ and[5Ji, where we show the radial profiles of the 
magnitude of the fluctuation (rms) value of the infall ve- 
locity together with the rms of the lateral component of 
the velocity and the rms value of the total velocity, re- 
spectively. We show in Fig. 2e, the radial profile of the 
Mach number which illustrates that the infalling veloc- 
ities are supersonic, while inside the core, the velocities 
are subsonic. 

Gravitational compression during the collapse of a pri- 
mordial gas cloud can at most lead to an amplification of 
the magnetic field strength by a factor of ~ p 2 / 3 in the 
limit of perfect flux freezing (i.e., ideal MHD). A stronger 
increase implies the presence of an additional amplifica- 
tion mechanism. Starting from an initial field strength 
of ~ InG, our simulations show a total magnetic field 
amplification by six orders of magnitude, leading to a 
field strength of about ~ 1 mG for the case where we 
resolve the local Jeans length by 128 cells. This is il- 
lustrated in Fig. [3K Fig. [3Jd shows that in our highest 
resolution simulation, the obtained field amplification is 
indeed stronger than what is expected from pure flux 
freezing, which demonstrates that the small-scale turbu- 
lent dynamo provides significant additional field ampli- 
fication over compression. Fig. [3J: shows a plot of the 
relative amplification in simulations with 16, 32, 64 and 
128 cells compared to the simulation where Aj is resolved 
by 8 cells, i.e., we divided each curve in Fig. [5Jd by the 
curve of the 8 cell run. The time evolution of the mean 
density p m within the central Jeans volume is depicted 
in Fig. [5H, while Fig. [5^ shows the plot of the rms ve- 
locity. The presence of turbulence delays the collapse 
until t ~ 4. During this time, the mean density shows 
oscillations while the rms velocity decreases as the tur- 
bulence decays. Sufficient numerical resolution is a cru- 
cial issue when studying the small-scale turbulent dy- 
namo. This is evident from Fig.Hl where the simulations 
with Aj resolved by 8, 16 and 32 cells are evolved fur- 
ther till t = 16. The 8 and 16 cell runs do not show 
any increase in B rms /p^l 3 and only exhibit weak fluc- 
tuations, while the 32 cell run is the first to show an 
increase. This signifies that a minimum of 32 cells per 
Jeans length is required to obtain the exponential ampli- 
fication of the magnetic field by the small-scale dynamo. 
Th is is also consistent w ith the turbulence simulations 
by iFederrath et ail (|2010f ). who concluded that at least 
30 grid cells are required to resolve turbulent vortices. 

As mentioned in Section 2, the growth rate of the 
magnetic field depends on the Reynolds number of the 
system and is thus related to the numerical resolution 
of the simulation. With increasing Reynolds number 
and thus with higher numerical resolution, the growth 
rate of the dynamo-generated magnetic field increases 
(|Haugen et alJl2004t iBalsara et all 12004). In agreement 
with this, our resolution study shows a divergent be- 
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Fig. 3. — Evolution of the dynamical quantities in the central 
Jeans volume as a function of r, defined in equation |(2]l for five 
runs with different number of cells to resolve the local Jeans length. 
Panel (a) shows the rms magnetic field strength B rms , amplified 
to 1 mG from an initial field strength of f nG, (b) the evolution 

of -Brms/Pm 3 j showing the turbulent dynamo amplification by di- 
viding out the maximum possible amplification due to perfect flux 

freezing, (c) the relative amplification in B rms /pm 3 compared to 
the 8 cell run, (d) the evolution of the mean density p m and (e) 
the rms velocity » rmB . The runaway collapse commences at about 
r ~ 4. 

havior in the growth rate of the magnetic field (Fig. \3]p 
and [3J;) rather than a convergent one. A convergence 
in the growth rate can only be obtained if the physical 
length scales that determine the Reynolds number are 
sufficiently resolved. Our simulations also show no signs 
of saturation and thus we simply stop the calculation 
when the numerical cost becomes prohibitively high. In 
reality, we expect the field amplification to continue until 
back- reactions either via the Lo rentz force (|Subramanianl 
U999tlSchekochihin et aT1l2004H o r via non-ideal MHD ef- 
fects such as ambipolar diffusion (|Pinto fc Gallil20"08l) be- 
come important. Calculations of MHD turbulence with- 
out self-gravity indicate maximum field strengths within 
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Fig. 4. — The figure illustrates the minimum resolution crite- 
rion required to capture the growth of the magnetic field due to 
small-scale dynamo action. The dynamo begins to be observed for 
simulations where Aj is resolved by a minimum of 32 cells. Simula- 
tions performed with the Jeans length resolved by either 8 cells or 
16 cells are decaying in nature with weak fluctuations. The verti- 
cal line indicates the values of B rms /Pm 3 obtained in the different 
resolution runs at r = 12.2. 

5% of the equipartition value (|Subramanianl [l999'). We 
note that the physical dissipation scales for ambipolar 
diffusion and Ohmic dissipation are much smaller than 
the Jeans length and thus the growth rates obtained in 
our simulations are lower limits to the physical growth 
rates. 

4. DISCUSSION AND CONCLUSIONS 

Studies of primordial star formation show that 

the gas within the first star forming ha los is 

turbulent (lAbel et al.1 [200l lO'Shea fc Normanl [20071: 
lYoshida et al.ll2008HTurk et al.ll2009f ). In this Letter, we 
showed that this turbulence is maintained throughout 
the collapse and drives a small-scale dynamo, which ex- 
ponentially amplifies weak magnetic seed fields. We per- 
formed five numerical simulations of collapsing Bonnor- 
Ebert spheres where the local Jeans length was resolved 
by 8, 16, 32, 64 and 128 cells. Our numerical simulations 
show that, starting from a weak seed field of ~ 1 nG, 
strong and dynamically important fields can eventually 
be generated in the central collapsing core during the for- 
mation of the first stars. The small-scale dynamo only 
works in simulations in which the turbulent motions in 
the central core are sufficiently resolved. We find that a 
minimum resolution of 32 cells per Jeans length is nec- 
essary for resolving the small-scale dynamo in our simu- 
lations. 

The generation of strong and dynamically important 
magnetic fields has interesting consequences for our un- 
derstanding of how the first stars form and how they 
influence subsequent cosmic evolution. We know from 
modeling galactic star-forming clouds that the presence 
of magnetic fields can reduce the level of fragmentation, 
and b y doing so strongly influences t he stellar mass spec- 
trum (jHennebelle fc Tevssierl l2008ft . Further more, the 
dynamo-generated strong magnetic fields can drive jets 
and outflows from the accretion disks via the ma gneto- 
centrifugal mechanism ([von Rekowski et~aTll2003[ ). Such 
outflows remove a significant fraction of the mass and 
angular momentum which influences the stellar mass 
spectrum. The implications of magneti c fields in self- 
gravitating disks have been explored by iFromang et al.l 
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(2004), who find that the interaction of turbulence ex- 
cited by the magneto-rotational instability with the self- 
gravitational instability excites additional modes, broad- 
ens the spiral arms and effectively reduces the accretion 
rate due to non-linear interaction. There are first at- 
tempts to stud y magnetic fields in the context of first 
star formation ([Machida et al.l [2006), but more sophisti- 
cated initial conditions and more appropriate magnetic 
field geometries need to be considered. 

Once the first stars have formed, they are likely to pro- 
duce a copious amount of ionizing photons, which drive 
huge Hn regions, bubbles of ionized gas, expanding into 
the low-density gas between the halos. The effects of 
radia tion have been shown in contemporary star form a- 
tion ([Krumholz et al.ll2005l [20071: iPeters etaLll2010allbh . 
and are likely to be important for th e first stars as well 
(I Wise fc Abell[2008tlGreif et al.l[200l . The expansion of 
the Hn region could be substantially different if mag- 
netic outflows drive a cavity into the surrounding gas. 
The magnetic field may further affect fluid instabilities 
near the ionization front. 

The mechanism of exponentially amplifying weak seed 
magnetic fields via the small-scale dynamo is likely to 
work not only during the formation of the first stars, but 
in all types of gravitationally bound, turbulent objects. 
Highly magnetized gas is thus expected already in the 
first galaxies. 
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